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current-voltage (I-V) characteristics, it is found that up to some constants I-V curves at a given temperature are 
identical to each other in a broad range of external currents. Simulations of the Villain model and the modified 
2D XY model allowing stronger thermal vortex fluctuations are also performed with RSJ type of dynamics. 
The finite-size scaling suggested in Medvedyeva et al. [Phys. Rev. B (in press)] is confirmed for all dynamic 
models used, implying that this finite-size scaling behaviors in the vicinity of the Kosterlitz-Thouless transition 
are quite robust. 
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I. INTRODUCTION 



The phase transition between the superconducting and 
normal states in many two-dimensional (2D) systems is of 
Kosterlitz-Thouless (KT) type The thermally excited vor- 
tices in a large enough sample, interacting via a logarithmic 
potential, are bound in neutral pairs below the KT transition 
temperature Tkt and as the temperature T is increased 
across Tkt from below these pairs start to unbind. The KT 
transition has been observed in experiments on superconduct- 
ing films 2D Josephson junction arrays and cuprate 
superconductors [gj. In these experiments, the current- voltage 
(I-V) characteristics have been commonly measured to de- 
tect the transition. For small enough currents it has a power 
law form V oc I a (or equivalently, E oc J a with the electric 
field E and the current density J), where the I-V exponent 
a is known to have a universal value 3 precisely at the tran- 
sition [Q]; For T < Tkt one has a > 3, whereas a = 1 
for T > Tkt [01- The dynamic critical exponent z which re- 
lates the relaxation time t to the vortex correlation length £ via 
r ~ £ 2 , is connected with a through the relation a = 1 +z [§]. 
Consequently, z has the value 2 at Tkt- 

Indeed, the values z = 2 and a = 3 at the KT tran- 
sition have been confirmed in many numerical simulations, 
e.g., the lattice Coulomb gas with Monte Carlo dynamics ||], 
the Langevin-type molecular dynamics of Coulomb gas par- 
ticles [|0p, and the 2D XY model both with the resistively 
shunted junction dynamics (RSJD) and the relaxational dy- 
namics (RD) ||ll],12| ]. That is why z » 6 obtained by Pierson 
et al. in Ref. [|D | at the resistive transition in many 2D sys- 
tems is very intriguing. We in Ref. Jl4] ] have re-analyzed the 
experimental I-V characteristics for an ultra thin YBCO sam- 
ple in Ref. [[jj and compared with the 2D RSJD model. This 
led to the suggestion that a novel finite-size type scaling ef- 
fect of the I-V characteristics (which can possibly be caused 
not only by the actual finite size of the sample but also by a 
finite perpendicular penetration depth or a residual weak mag- 
netic field [@Jl4p) is responsible for the large value of z w 6 



ponent which at Tkt has value two [14|. Since similar con- 
clusions were also reached in Ref. [ 15 1, we believe that there 



now seems to emerge some consensus [ J14UL5I ], although the 
physical mechanism causing this finite-size scaling behavior 
is still unclear. 

In the present paper we search for a possible reason of this 
finite-size scaling behavior. We use the RSJD [ pTJl^|i"oj ], the 
RD (often referred to as time-dependent Ginzburg-Landau dy- 
namics) [|l l[|l2|], and the Monte Carlo dynamics (MCD) sim- 
ulations |]9|, |17[|18| ] to study the I-V characteristics of three dif- 
ferent models of 2D superconductors: the Villain model fll^], 
the XY model in its original form, and the XY model modi- 
fied with p-type of potential (see Ref. [ ^0| ] for details). These 
three models differ from each other by the density of the ther- 
mally created vortices. We come to the conclusion that the 
qualitative features of the scaling suggested in Ref. Jl4] ] are 
ubiquitous and independent of both the vortex density and the 
type of dynamics. We also demonstrate that the I-V curves 
obtained with different types of dynamics, up to constant scale 
factors, coincide very well in a broad range of the external cur- 
rent density. This opens a possibility to use MCD simulation, 
which from a simulation point of view is more efficient than 
RSJD or RD, to examine long-time dynamic properties of the 
models. 

The layout of the paper is as follows: In Sec. || we reca- 
pitulate the Hamiltonians for the generalized 2D XY model 
with the p-type of potential (the usual XY model is recovered 
when p = 1) and the Villain model, and describe details of 
the dynamics used (RSJD, RD, and MCD). The results from 
the usual XY model with p = 1 subject to different dynamics 
(RSJD, RD, and MCD) are presented in Sec. |l| while the re- 
sults from RSJD simulations applied to the different types of 
models, i.e., the usual XY model with p = 1, the XY model 
with p = 2 and the Villain model, are described and analyzed 
in Sec. |[v|. We summarize and make final remarks in Sec. [v| 



II. MODELS AND DETAILS OF DYNAMICS 



in Ref. [Ofl, but at the same time it was concluded that this 
exponent has nothing to do with the true dynamic critical ex- 
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A. Model Hamiltonians 

The 2D XY model defined on a square lattice, where each 
lattice point i is associated with the phase 9i of the supercon- 
ducting order parameter, is often used for studies of the KT 
transition. The phase variables in this model interact via the 
Hamiltonian with nearest neighbor coupling, which in the ab- 
sence of frustration is given by 



H = J2 Ufa = Oi - 9j), 



(1) 



where (ij) denotes sum over nearest neighbor pairs, <fi is the 
angular difference between nearest neighbors, and the inter- 
action potential U (4>) is written as 



17(0) = Ej(l 



(2) 



with the Josephson coupling strength Ej. 

The dominant characteristic physical features close to the 
KT transition are associated with vortex pair fluctuations. One 
interesting aspect is then how the density of the thermally ex- 
cited vortex fluctuation effects the critical properties. To study 
this we generalize the interaction potential by using a param- 
eter p pOj] : 



U*{<j>) = 2Ej 



1 — cos 



2,r 



(3) 



where U p=1 (4>) corresponds to the potential of the usual XY 
model [see Eq. (Q)]. The practical point with such gener- 
alization is that the vortex density increases with increasing 
p J2C|]. The variation of the parameter p can also change the 
nature of the transition: for p exceeding some maximum value 
(p > Pmax ~ 5) the t ype of the phase transition changes from 
KT to the first order [20,21]]. In the present paper we choose 
p = 2 which is well inside the KT transition region, yet is 
large enough to ensure substantially more vortex fluctuation 
over a temperature region around the phase transition in com- 
parison with the usual p — 1 XY model given by Eq. (^). 
While the XY model with the p-type potential with p > 1 
has more vortices than the usual p = 1 XY model, we also 
study the Villain model jl^ ] which has less vortex-antivortex 
pairs [p2[|. The interaction potential U(<p) in the Villain model 
is given by 



,-UW/T 



exp 



E 



2T 



J -{2m-(j)f 



(4) 



To simulate the dynamic behaviors of these models we use 
several types of dynamics: RSJD, RD, and MCD. All these 
dynamics should result in the same equilibrium static behav- 
iors if we apply them to the models with the same interaction 
potential. However, dynamic properties of the systems can be 
different. Of course, different types of dynamics have their 
own advantages and disadvantages. The RSJD is constructed 
from the elementary Josephson relations for single Josephson 
junction that forms the array units, plus Kirchhoff's current 



conservation condition at each lattice site | |1 1|J12| ]. Therefore, 
this type of dynamics has a firm physical realization. On the 
other hand, RSJD is quite slow which leads to the limitation 
in the time scale one can probe in simulations. Although the 
RD JTT| , p^pT1 | is much easier to implement than RSJD, it does 
not converge much faster than RSJD and it does not have a 
similar direct physical realization as RSJD. However, a su- 
perconductor has been argued to have a RD type of dynam- 
ics rather than a RSJD [gQ. The MCD simulations ^ 
are much faster than RSJD or RD, which allows one to in- 
vestigate dynamic behaviors in much longer time scale (one 
can also study dynamic behaviors at much lower temperatures 
with MCD). However, since there is no direct physical real- 
ization of the MCD in practice the applicability of this dy- 
namics to a specific physical system must then be explicitly 
demonstrated. In the following discussions on the details of 
the different dynamics used, we focus on the original 2D XY 
model with the interaction potential in Eq. (||) since the exten- 
sions to a modified 2D XY model (||) and Villain model (Q) 
are straightforward. 



B. Dynamic models 



In this section we briefly review the dynamical equations 
of motion for RSJD, RD, and MCD, in the presence of the 
fluctuating twist boundary condition (FTBC) [[H]|l8|]. We 
perform simulations of unfrustrated square L x L lattices 
with L = 6, 8, and 10 at various temperatures to mea- 
sure the voltage across the lattice as a function of the exter- 
nal current. Although the system sizes are relatively small, 
which is inevitable because of the low temperatures and 
the small external currents used here, the FTBC has been 
shown to be very efficient in reducing the artifact due to 
small system sizes JTT| , p4| ], and reliable results can be estab- 
lished [yJHH]. 

In the FTBC, the twist variable A = (A x ,A y ) is in- 
troduced and the phase difference ^ij on the bond is 
changed into B. L — 9j — • A, with the unit vector 
from site i to site j, while the periodicity on L is imposed: 
9% = &i+L± = ®i+Ly- The Hamiltonian of 2D L x L XY 
model under FTBC without external current has been intro- 
duced in Ref. [0], and is written as [compare with Eq. ([[])] 



(ij) 



cos 



- 6j - r tf • A), 



(5) 



which later in Ref. [ |1 8| ] has been extended to the system in the 
presence of an external current and written as 



H = -Ej c °s &j + ^-L 2 JA X , 

do) 



(6) 



where J the current density in the x direction. 
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1. RSJD and RD 

We introduce first the RSJD equations of motion for phase 
variables and twist variables, which are generated from the 
local (global) current conservation for the phase (twist) vari- 
ables (see Ref. [nil] for details and discussions). The net cur- 



rent I*.-. Iij 



rent Iij from site i to site j is the sum of the supercurrent Ifj , 

the normal resistive current and the thermal noise cur- 

l j 

I?j + ^ij + tfj' The su P ercu iT en t is given 
by the Josephson current-phase relation, Ifj = J c sin<^y, 
where I c = 2eEj/h is the critical current of the single junc- 
tion. The normal resistive current is given by = Vij /r, 
where Vij is the potential difference across the junction, and 
r is the shunt resistance. Finally the thermal noise current 
I\j in the shunt at temperature T satisfies (ij.) = and 
<4(^«(0)) = {2k B T/r)8{t){8 ik 8ji - 6 u 6 jk ), where (■ • ■) 
is thermal average, and S(t) and Sij are Dirac and Kronecker 
delta, respectively. Using the current conservation law at 
each site of the lattice together with the Josephson relation 
4>ij = d<fiij/dt = 2eVij/h one can derive the RSJD equations 
of motion for phase variables: 



(7) 



where the primed summation is over the four nearest neigh- 
bors of j, Gij is the lattice Green function for 2D square lat- 
tice, and rjjk is the dimensionless thermal noise current de- 
fined by rjjf. = l!j k /I c . The time, the current, the distance, 
the energy, and the temperature are normalized in units of 
Ti/2erl c , I c , the lattice spacing a, the Josephson coupling 
strength Ej, and Ej /&b, respectively. In order to get a closed 
set of equations we further specify the dynamics of the twist 
variable A from the condition of the global current conserva- 
tion that the summation of the all currents through the system 
in each direction should vanish [nil : 



L 2 ^ S1I1( W 



A v = J2 Sin< 



VA W - J, 



(8) 
(9) 



where Yl/ij) denotes the summation over all nearest neigh- 
bor links in the x direction, and we apply the external dc 
current with the current density J in the x direction. Here, 
the thermal noise terms r]A a and rjA v obey the conditions 
(i7A„) = (»7A B ) = (VA x VA y ) = 0, and (r)A x (t)VA x {°)) = 
( V A y (t)r lAy (Q))^(2T/L^5(t). 

In the RD, the equations of motion for the phase variables 
are written as [fll|l2LE6|| 



% = -r 



.OH 
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sm < 



(10) 



where T is a dimensionless constant (we set T = 1 from 
now one), H is in Eq. (0), t is in units of h/TEj, and 



the thermal noise rji(t) at site i satisfies (r]i(t)) = and 
('7*O0 ? 7j(O)) — 2T5(t)5ij. The equation of motion for the 
twist variables in the absence of an external current is of the 
form (see Ref. [|ll|,|l2|] for more details) 



1 dH 



?7A, 



(11) 



which is the same as Eqs. (|8j) and (Q) for RSJD. Accordingly, 
to some extent the RD may be viewed as a simplified version 
of the RSJD where the global current conservation is kept but 
the local current conservation is relaxed. 

Consequently, the equations for the phase variables are dif- 
ferent for RSJD and RD [Eqs. (Q) and (|l0|), respectively] while 
the same equations (||) and (^) apply to the twist variables for 
both dynamics. These coupled equations of motion are dis- 
cretized in time with the time step A< = 0.05 and At = 0.01 
for RSJD and RD, respectively, and numerically integrated 
using the second order Runge-Kutta-Helfand-Greenside algo- 
rithm [^7p. The voltage drop V across the system in the x 
direction is written as V = —LA X (see Ref. jTl]]) in units 
of I c r for RSJD and in units of TEj/2e for RD, respectively. 
We measure the electric field E = (V/L) t to obtain I-V char- 
acteristics, where (•••}* denotes the time average performed 
over O(10 6 ) time steps for large currents for both RSJD and 
RD, and 0(5 x 10 7 ) and O(10 9 ) steps for small currents for 
RSJD and RD, respectively. 



2. MCD 

The technique to simulate 2D XT model with MCD is based 
on the Hamiltonian (^) and the standard Metropolis algo- 
rithm [[?8]]. The one MC step, which we identify as a time 
unit, is composed as follows [18]]: 

1 . Pick one lattice site and try to rotate the phase angle at 
the site by an amount randomly chosen in [—66, 56] (we 
call 56 the trial angle range). The twist variable A is 
kept constant during the update of the phase variables. 

2. Compute the energy difference AH before and after the 
above try; If AH < or if e ~ Aff / T is greater than a 
random number chosen on the interval [0, 1), accept the 
trial move. 

3. Repeat steps [TJand|2]for all the lattice sites to update the 
phase variables. 

4. Update the fluctuating twist variables A x in the similar 
way that LA X is tried to rotate within the angle range 
5 A with 6i and A y kept unchanged. (For convenience, 
we use 5 A = 56). 



5. Compute the energy difference AH before and after the 
trial step [| for 
erwise accept it with ] 



y umeiciiuc L-±n uciuie anu anci inc 

Accept the step EL if AH < oth- 
h probability e -^™/ T like in step |[ 



6. Repeat steps EJ and || to update A y . 
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In the MCD simulation the trial angle 69 = tt/6 has been 
chosen since it is sufficiently small in order to obtain the cor- 
rect I-V characteristics while it is big enough to make MCD 
much faster than the other dynamic methods [|l8]| . The time- 
averaged electric field is obtained after equilibration from the 
averages over O(10 9 ) (at large currents) to O(10 10 ) (at small 
currents) MC steps. 



III. 2D XY MODEL SUBJECT TO DIFFERENT TYPES OF 
DYNAMICS 

In this section we use three different types of dynamics, the 
RSJD, the RD, and the MCD to study dynamic behavior of 
the 2D XY model with p = 1 under the FTBC. The dynamic 
behavior of the system can be obtained from the complex con- 
ductivity, the flux noise spectrum, as well as the I-V charac- 
teristics which is commonly measured in experiments. We 
will focus on the I-V characteristics in the present paper. As 
pointed out in Sec. |l| the RD to some extent may be consid- 
ered as a simplified version of the RSJD. Thus from this point 
of view it is perhaps not surprising that these two models (as 
we will see) contain similar features of the vortex dynamics. 
In Ref. [^] from the study of a simple dynamic model of iso- 
lated magnetic particles in a uniform field it has been shown 
that the actual dynamics of the model and MCD are in a good 
agreement when the acceptance ratio of the Metropolis step is 
low enough. This implies that the MCD should give the same 
I-V characteristics as the RSJD after an appropriately chosen 
normalization of time, when the trial angle 69 is sufficiently 
small (it was shown in Ref. [ 18 1 that 69 = ir/6 is sufficiently 



small). We will confirm this further in the present simulations. 

In Fig. [I] we compare I-V characteristics in the form of 
the electric field E versus the current density J obtained from 
RSJD, RD, and MCD simulations in the temperature range 
0.70 < T < 1.50 (the temperature interval is 0.05 if 0.70 < 
T < 1.00 and 0. 10 otherwise) in log scales for the system size 
L = 8. In order to make I-V data of RSJD and RD simula- 
tions coincide, E obtained with RD (Erd) is multiplied by a 
temperature-independent factor represented by the horizontal 
line in the inset of Fig. [I]. Since the measured time-averaged 
electric field is inversely proportional to the time scale for a 
given dynamics, one can from the ratio Ersjd /Erd infer 
the correspondence between times of RD and RSJD. ^From 
this comparison, we find for 0.70 < T < 1.00 that one unit 
of time in RD approximately corresponds to 0.526 time unit 
in RSJD, independent of the temperature. Note that the I- 
V curves corresponding to the temperatures exceeding 1.00 
are almost a straight lines. Therefore the collapse of data be- 
tween the different dynamics is trivial for T > 1.00. Also 
the I-V characteristics obtained from MCD can be made to 
collapse on top of the corresponding curves for the RSJD and 
RD, as shown in Fig. [l] However, in this case the time scale 
factor, describing how many RSJD time steps one MCD step 
corresponds to, depends on the temperature, as shown in the 
inset of Fig. [j], where this factor is shown to be a linear func- 
tion of temperature for system size L = 6, 8, and 10 up to 
T = 1.00. This is in accordance with the model studied in 



Ref. [|29p, where the same linear behavior in terms of the tem- 
perature has been found. Thus for the I-V curves we have 
a precise relation between RSJD and MCD: For example, at 
T = 0.80 we get 1 MCS = 10.2xRSJD time unit. This opens 
a practical possibility to study some aspects of dynamic be- 
havior of the XY model by using MCD simulation, which is 
usually much more efficient than RSJD and RD. 

We have shown in this section that the RSJD, the RD, 
and the MCD applied to the usual XY model gives basically 
the identical I-V characteristics up to some constant factors. 
From this observation, one can conclude that the dynamic crit- 
ical behaviors of the XY model inferred from the I-V charac- 
teristics should be identical for all these dynamic models. We 
in next section use the RSJD to study the Villain model and 
XY models with p — 2 and presume from the observation in 
the present section that the conclusion drawn in Sec. [v] for 
the RSJD case should be also valid in the other dynamics (RD 
and MCD). 



IV. VARIOUS XY MODELS SUBJECT TO RSJD 

To study the critical behavior of the system in the vicinity 
of the transition one can use scaling relations. Fisher, Fisher, 
and Huse (FFH) in Ref. proposed that the nonlinear I-V 
characteristics in a D-dimensional superconductor scales as 



E = J£ 



D-2-z 



(12) 



where £ and z are the correlation length and the dynamic crit- 
ical exponent, respectively, and \± is the scaling function 
above (+) and below (— ) the transition. Pierson et al. in 
Ref. [|l3|] have applied a variant of this FFH scaling approach 
to the I-V data for thin (D = 2) superconductors and su- 
perfluids and suggested a phase transition with z w 6. In 
Ref. [14 1 another scaling relation has been introduced and it 
has been shown that a certain finite-size effect which is not in- 
cluded in the FFH scaling may have caused the large z in spite 
of the fact that the finite size effect precludes the possibility 
of a real phase transition. This finite-size scaling around and 
below the KT transition is given by the form (see Ref. [ |l4| ] for 
the details) 



E 

Jr 



LlLg L (T)j, (13) 



where R — \imj^ (E/J) oc L^ z ^ is a finite-size induced 
resistance without external current, /it(0) = 1, h(x) oc x z ^ 
for large x, and ffi(T) is a function of at most T and L 
such that a finite limit function goo(T) exists in the large-L 
limit. For small values of the variable x = JLgjj(T) the In- 
dependence of the scaling is absorbed in a function g^ (T) for 
each fixed size L giving rise to the scaling form 



^ = h(jLg L (T)). 



(14) 



In Fig. ^ we demonstrate the existence of the finite-size 
scaling given by Eq. ( |l4| ) for L = 8 within the temperature 



4 



intervals 0.70 < T < 1.10. The data are obtained for MCD 
and scaled by the appropriate factor so as to correspond to 
RSJD and RD (compare Fig. |lj). Because of the correspon- 
dence between the three types of dynamics (see Sec. this 
also means that the existence of the finite-size scaling given 
by Eq. ( |l4] ) is insensitive to the choice of dynamics. 

The FFH scaling given by Eq. ( |l2| ) is correct only in the 
thermodynamic limit £/L — » 0. However, from a practical 
point of view there is a connection between Pierson method, 
which is based on FFH scaling, and the finite-size scaling 
introduced by Eq. (|l4|). If we assume that £ in Eq. ( |l2] ) is 
proportional to R~ a , where a is T-independent constant, the 
connection between the two different scaling approach be- 
comes: 3l(T) = Aj,R~ a with Al being a constant which 
may depend on L. In Fig. || there are presented three different 
functions gh{T) corresponding to L = 6, 8, and 10. These 
functions are determined from the condition that curves cor- 
responding to the different temperatures should collapse when 
plotted as E/ JRvs LJgL(T). Since it is well established that 
the 2D XY model on the square lattice has the KT transition 
at T KT w 0.892 (Ref. [0), Fig. | shows that g L {T) over a 
limited region in the vicinity of the KT transition is very well 
represented by the R~ a with a = 1/6 for all investigated 
system sizes. 

Next we demonstrate how the finite-size scaling given by 
Eq. ( p4[ ) works for the I-V data obtained by simulations of the 
modified XY model. These simulations are done with RSJD. 
Fig. ^(a) verifies that this scaling indeed exists for the XY 
model modified with the Villain type of potential introduced 
by Eq. (|]). The inset of Fig. Q(a) shows the scaling function 
gh{T) determined by finding the best data collapse for small 
values of LJg^(T). One can see that in the vicinity of Tkt> 
which is approximately equals to 1.4 for this model, 3l(T) 
can be fitted by Ar,R~ a with a = 1/6. The data collapse in 
Fig. [|(b) shows that E/ JR is only a function of the scaling 
variable JLgL(T) (when JLgL{T) is small enough) for the 
XY model modified with p-type of potential (Eq. (Q)), where 
p = 2. The inset shows the function <?l(T) together with 
AlR~ u - Since Tkt ~ 1.15 for this model one can again see 
that the scaling function <7l(T) in the vicinity of KT transition 
is proportional to R~ a with the same exponent a = 1/6 as for 
the original XY and Villain models. 

The crucial difference between the Villain model, the usual 
XY model and the p = 2 XY model in the present con- 
text is the vortex density. The KT-transitions for these models 
occur at the Coulomb gas temperatures T GG = 0.23, 0.2, 
and 0.1, respectively (T CG = ^(U") (see Ref. g)). Lower 
T GG means higher vortex density. Thus the finite-size scal- 
ing property given by Eq. ( |l4| ) appears to be independent of 
vortex density. 



V. DISCUSSION AND CONCLUSIONS 

We have simulated 2D XY model with three types of dy- 
namics: RSJD, RD, and MCD. The main conclusion of the 
paper is that the qualitative features of the finite-size scaling 



given by Eq. ( |l4"| ) are independent of both the vortex density 
and the type of dynamics. Therefore the finite-size scaling be- 
havior given by Eq. @ of the finite-size induced tails of the 
I-V characteristics appears to be a robust feature. From the 
comparisons of the current-voltage characteristics obtained 
for each type of dynamics we found that, up to some scale 
factor, I-V curves at a given temperature are identical over 
a broad range of external currents. This makes it possible to 
use MCD simulations, which are more computer efficient than 
RSJD and RD simulations, to obtain I-V curves correspond- 
ing to both RSJD and RD dynamics. 

The phase transition for the 2D XF-type models are of KT 
type with z = 2. This raises the intriguing question of the 
origin of the large z w 6 obtained by Pierson et al. [|l3|]. In 
Ref. [ |l4] ] it was argued that the Pierson scaling in relation to 
the finite-size scaling given by Eq. ( |l4] ) corresponds to the 
proportionality gL(T) ~ R a where 1/a is the exponent 
which corresponds to the "z" obtained by the Pierson scal- 
ing. The reason for the existence of this scaling like behavior 
is still unclear. 

In the present paper we have shown that, within the class of 
2D XY-type models studied, a value 1/a w 6 is obtained in- 
dependently of the type of dynamics, as well as, of the vortex 
density. The origin of this seemingly robust behavior calls for 
further investigations. 
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FIG. 1. Comparisons of the I-V characteristics for the 2D XY 
model with MC, RSJ, and RD dynamics on a square lattice with the 
finite size L — 8 plotted as E — V/L against J = I/L in log 
scales at temperatures T = 1.5, 1.3, 1.2, 1.1, 1.0, 0.95, 0.90, 0.85, 
0.80, 0.75, and 0.70 (from top to bottom). Multiplication by a fac- 
tor (presented in the inset), which in the case of MCD depends on 
temperature, makes the curves coincide in a broad range of external 
current density. The straight line from the origin in the inset shows 
that the factor in the MCD case is a linear function of temperature 
for system sizes, L — 6, 8, 10 in accordance with expectation |p9||. 
However, for higher T there is a deviation. The horizontal full line in 
the inset shows that the factor between the RSJD and RD I-V curves 
is independent of T (the factor is multiplied by 10 in the inset). 
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0.01 0.1 1 10 100 

Ug L (T) 

FIG. 2. The finite-size scaling of the current-voltage characteris- 
tics, E/JR vs LJg L (T), given by Eq. (Q. The function g L (T) 
is determined from data collapse for small values of LJgh(T). The 
data was obtained from MCD and was scaled with a T-dependent 
factor so as to correspond to RSJD and RD; system size L = 8 and 
1.10 < T < 0.70 where each curve corresponds to a fixed T. 



FIG. 3. The function <?i(T) determined for I-V characteristics 
with size L — 6 (open squares), L = 8 (open circles), L = 10 
(open triangles) obtained as in Fig. Q One notes that the function 
gL for these sizes over a limited T interval around KT-transition (at 
T ~ 0.90) is well approximated by gL oc R~ a with a ~ 1/6. The 
data for L = 10 function <?l(T) as well as AlR^ 01 are multiplied 
by 1.5 for convenience. 
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FIG. 4. Existence of the scaling in the form Eq. ( |14| ) for modified 
2D XY models corresponding to different vortex density. Systems 
with size L — 8 have been simulated with RSJD. The finite-size scal- 
ing of the I-V characteristics obtained for 2D XY model with the 
Villain potential in the temperature range 1.70 < T < 1.15 is shown 
in (a) whereas (b) shows the same thing for the p = 2-potential for 
1.30 < T < 0.95. The inset on both (a) and (b) shows the function 
gi(T) determined from the condition of the best data collapse. It is 
also shown there that the function gL for both cases over a limited T 
interval is well approximated by gL oc R~ a with a « 1/6. In (a) 
g L oc 0.76i? _1/6 and in (b) g L oc 1.04ir 1/6 . 
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